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Abstract. Optical photometry 
PSR0531+21 has been extended 
tions that strengthen evidence for a 
60 seconds periodicity. This period 
creasing with time at approximately 
rotational period of the pulsar. The 
its time dependence fit a simple free 



for the pulsar 
with new observa- 
previously observed 

is found to be in- 
the same rate as the 
observed period and 
precession model. 



Key words: Stars: neutron - pulsars: individual: 
PSR0531+21 - rotation - oscillations 



1. Introduction 

More than seven years ago we began an investigation of 
the optical light curve of the Crab pulsar. We tentatively 
identified a peak in Fourier spectra of the light curve at 
a freque ncy « 1/60 Hz (Cadez & Galicic 113561 a, IT5561 
b, 1996 c). Later observations slightly enhanced the signal 
to noise ratio of the suggested peak. We suspected that it 
was the signature of free precession of the pulsar and thus 
in EEO proposed a simple model. We tested it by predict- 
ing the free precession frequency of the Earth and showed 
that it also reproduces well the measured ellipticities of 
the planets (Cadez et al !1997|l . The former test was based 
on the argument that the relative rigidity of the pulsar is 
at least as strong as that of the Earth. This allowed us 
to apply the same model to calculate the free precession 
frequency of pulsars. We found that the observed 1/60 Hz 
frequency is consistent with a 1.3 Mq pulsar model based 
on the tensor interaction equation of state. As a further 
check we proposed to measure the slow-down rate of free 
precession. In the past years we gathered more data, build- 
ing a set of more than 3.400 images, which cover more than 
20 hours of photometry spanning a period of almost nine 
years. It is appropriate to confront these data with the 
proposed test. 



2. Calibrations and data sets 

Our photometry was done in the stroboscopic mode. This 
type of photometry was introduced previously (Cadez & 
Galicic 1996) in order to increase the signal to noise ratio 
of measured pulsed emission from pulsars. The method is 
based on the fact that the arrival time and duration of 
pulses can be calculated quite precisely. Thus we intro- 
duced the stroboscope - a phase controlled shutter (rotat- 
ing wheel), which passes light to the CCD detector only 
during the time intervals when the pulsar (main pulse) is 
"on" . Therefore, images obtained behind this shutter re- 
ceive all light emitted by the pulsar, but the much brighter 
emission from the surrounding nebula, is reduced by the 
ratio of the expected main pulse duration to the pulse pe- 
riod (set to 0.1 by our stroboscope). In our case this noise 
reduction ratio is about a factor of 2. In stroboscopic pho- 
tometry the timing is crucial, so we spent a considerable 
effort checking it. The whole stroboscopic setup together 
with the CCD has been periodically tested at the Hewlett- 
Packard calibration laboratory in Ljubljana, where an ar- 
tificial pulsar has been pulsed by a cesium clock (Galicic 
1999}. All tests confirmed that the stroboscopic wheel cor- 
rectly follows the frequency and phase of the input signal 1 . 

Photometric data were taken with four different tele- 
scopes and exposure times between 4 and 15 seconds 2 and 
sampling rates beween 4.0 and 26.7 seconds. The data are 
of a somewhat varying quality depending on the telescope 
and on respective observing conditions. In Table 1 we list 
the dates, telescopes, durations of observation and sam- 
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1 During the observations abbreviated by As we had some 
timing problems which were traced to a software bug in cal- 
culating the proper date, so that the calculated barycentric 
frequency set on the timer was slightly off, resulting in slow 
phase slippage of the stroboscope with respect to the pulsar. 
The slippage was slow enough that we could calculate its effect 
and correct for it. 

2 The 4 seconds sampling rate of the HST data is not intrinsic 
to the high speed photometer but rather to our processing 
procedure. 
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Fig. 1. Light curves of the compari- 
son star. 



pling rates of all runs. The abbreviations HST (Percival 
et al. I1993[) CE, As and GH stand for Hubble space tele- 
scope, the 1.82 m telescope of the Asiago observatory, the 
1.22 m telescope of the Asiago observatory and the 2.12m 
telescope of the Guillermo Haro observatory respectively. 

In Fig. 1 and 2 we show the light curves of a field star, 
which is used for calibration and the contemporeneous 
light curves of the Crab pulsar. The photometry was done 
by using IRAF/DAOPHOT and previously DAOPHOT 
routines 3 . 



3. Data analysis and results 

In order to search for a possible common modulation fre- 
quency in all Crab light curves, we construct a matrix 
with components Mjk(u>), which are Fourier transforms 
of cross-correlation functions between light curves, i.e.: 

M jk {uj) = \Fj(id)F k (uj)\ , (1) 
and Fk{to) is the Fourier transform of the k th light curve: 



considered as a continuous function defined first at dis- 



(2) 



Here to/, (ti ) is the magnitude as measured from the I im- 
age of the k th light curve taken at time ti (with respect to 
the beginning of the k th set) and T k = i; mQX — tj m4n is the 
total duration of the k th set of measurements. Fk(u>) is 



3 Until 1998 the DOS version of DAOPHOT available at 
www.fiz.uni-lj.si/astro/daophot.html and later IRAF V2.ll 
has been used. 



crete frequencies uj s 



where s is an integer be- 



longing to the interval < s < {l ma x — Imin)/^, and 
then extended so that |i*)-(u;)| = |Ffc(o; s )| on the interval 

We first test statistical properties of the Mjk compo- 
nents on light curves of the field star, i.e. we test the hy- 
pothesis that all noise in the light curve is white Gaussian 
noise. If this hypothesis is correct, then (for each ui) Mij 
(i 7^ j) may be considered as a realization of a random 
process /i distributed according to: 



"to -$*>($) 



(3) 



and Ma belongs to a random process v distributed ac- 
cording to 



(4) 



where 4 a 2 — (\Fk(uj)\ 2 ) is the spectral density of noise in 
measuring the magnitude, Kq(x) is the modified Bessel 
function (Galicic 1999) and () indicates the average over 
frequencies (higher than 0.005Hz to exclude the 1/f noise 



4 The average is taken over all k and uj and, since noise is 
considered as white, all Fk belong to the same random process. 
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at the low end). In order to detect a periodic component 
in the light curve we construct two quantities 5 : 



M(lo) 



and 



1 



N N 



i=i j=i 



(5) 



\ 



N N 



i=l j=l 

which in case of white Gaussian noise average to: 



(MM) = ^-i^V 



2N(N+l) 



(6) 



(7) 



and: 

(o-m(cj)) = 



(8) 



1/(16 - ir 2 )N 2 + 2(32 + 7T 2 - 8tt)A^ + (16tt - 16 - tt 2 ) 



2(JV+1) 



In Fig. 3 we plot M(uj), ctm(w) and ^f^rjj calcu- 
lated from the 12 light curves of the field star as plot- 
ted in Fig. 2. Here cr' M (cj) = <Jm{oj) if 0m(w) > (ctm(u)) 
and o-' m (uj) = (<t m {oj)) if cm(w) < (ffjf(w)}. The aver- 



age <M(w)) = 0.00091^ff and (ctmM) = 0.00054^- 



5 In fact, since the light curves are of varying quality, we 
define M and om as weighted averages, where the weight of 
each component Mjk(u) is -^t^t, where CTj(k) is the average of 

(cr(u))) corresponding to the j th (k th ) light curve. 



Table 1. List of data sets giving the name of the run 
(as used in later references), the name of the telescope, 
the date of observation, the sampling rate (the exposure 
time to obtain the sample is generally shorter than the 
sampling rate) and the total duration of the run. 
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Fig. 2. All pulsar light curves. The relative magnitude 
scale is at lower right. 



Thus, from eq|7|we deduce a = 0.023^^ 
a = 0.020^£ 

/Hz 



and eq|Sl gives 
Both numbers are consistent with each 
other and with average dispersion estimates from IRAF 
{o~iraf = o-/^/t; t = exposure time). Therefore, we con- 
clude that the test star data conform with the statistical 
assumptions described above. 

The result of a similar analysis for 16 Crab light curves 
is shown in Fig. 4. It is apparent that (ctm(w)) is almost 
the same for the Crab and the test star; however the av- 
erage (M(u)) for the Crab is 1.5 times that of the test 
star. The process Mjj(u) can, therefore, not be consid- 
ered as white Gaussian noise, but it can be understood 
by assuming that it is the sum of a stationary "signal" 
MfAw) and white Gaussian noise My (a;) with the same 
spectral density as that of the test star. We considered 
the possibility that the "signal" is connected to a spuri- 
ous phase modulation in stroboscope timing, but lab tests 
have set the upper limit for such a modulation some 100 
times lower (Galicic 1999). We therefore conclude that the 
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median ratio of the two (the upper 
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Fig. 4. The average cross- 
correlation function, its sigma 
and the median ratio for pulsar 
data. Symbols used are the same as 
those in Fig. 3. 



Crab pulsar light curve is intrinsically intensity modulated 
at the level of about 0.03 . The signal to noise ratio 

V Hz 

of this modulation measure is presently too low to reveal 
more about its nature. 
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Fig. 5. The average cross- 
correlation function (M(u))), its 
sigma ((Tjvf(o))) and the median 

ratio ( ^ ^ for pulsar data as 

the function of the recalibrated 
frequency v = £ = £ (^f), (i.e. 

p = 1). Other symbols are the same 
as in Fig 3. 

To test for the slow down of the purported free preces- 
sion we recalibrate frequencies in Fk{u>) to the same date 
according to an assumed power law dependence of the free 
precession frequency with time, i.e. we introduce the re- 
calibrated frequency u — u>(-^-) p , where v re f and Vk are 
the rotation frequencies of the pulsar at the reference date 
(Jan. 1.1996) and at the date when the k th light curve was 
obtained. The above analysis has been repeated for differ- 
ent powers p between and 4. A 5.3 a peak stands out in 
M (u))/a' M (u>) and is found at ui = s9 q 3s for p around 1, 
as shown in Fig. 5 and 6. Note that in the free precession 
model p is expected to be 1 if the time for the pulsar to 
relax its internal stress is much longer than the total ob- 
servation period (« 9years). It is expected to be 3, if this 
relaxation time is considerably shorter. Short relaxation 
times would point to externally driven precession. 

4. Discussion 

Thc question whether pulsars free precess was raised soon 
after pulsars were identified as rotating neutron stars. 
However, more than twenty years later most astronomers 
would agree with Trimble and McFadden's (1998) some- 
what humorous reference to our paper suggesting that 
Earth is the only known free precessing body in the Uni- 
verse. We believe that pulsar free precession has evaded 
detection and/or identification for two main reasons. 1) 
Most researchers assume that neutron star matter is in 
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a complex superfluid state, which provides an excuse to 
evade the question: what is the inertial eccentricity of a 
rotating neutron star? As a result, the range of predicted 
or "identified" free precession frequencies in the literature 
is almost without bounds, although Pines and Shaham 
l|1974(l estimated that the Crab pulsar would wobble with 
a period ~ 5 minutes if it had a solid core. 2) As Alpar 
and Ogelman i|1987[l correctly point out, free precession 
is a motion characteristic of rigid rotators. In a plastic 
rotator, free precessing motion is damped out by viscous 
losses - in pulsars this should happen on a time scale of 
years or thousands of years. Therefore, unless pulsars are 
exposed to reasonably large external torques, they are not 
expected to free precess. In our earlier paper (1997) we 
tried a naive approach: starting with the observation that 
the Earth free precesses but taking into account the pos- 
sibility that it may not be quite rigid with respect to free 
precession. We analysed the relation between the Earth's 
inertial eccentricity, free precession frequency, rotation fre- 
quency and obliquity. We found that the Earth's shape 
and that of other planets is quite similar to the equilibrium 
shape of rotating polytropes with no shear stress. The 
Earth's free precession frequency can be fairly accurately 
calculated on the basis of these assumptions. We argue 
that a pulsar's crust is relatively more rigid and thick than 
the Earth's. Therefore, argument for planets can also be 
used to calculate a pulsar's free precession frequency. For 
the fast Crab pulsar we found a free precession frequency 
on the order of minutes" 1 and in particular the theoret- 
ical value 59 l sec was found for a 1.3 M Q pulsar model 
based on the tensor interaction equation of state. Surpris- 
ingly, the 35 day period of HerX-1 fits the same formula 
almost exactly 6 . Our much larger dataset for the Crab 




Fig. 7. The spectral amplitude of 
the signal at the expected frequency 



pulsar strenghtens evidence for the 60 seconds period and 
also suggests that this period is increasing with time at 
almost the same rate as the rotation frequency of the pul- 
sar. We also found that different combinations of smaller 
data sets consistently produce a peak at lu = 5g 2 g 3s (of 
course with a correspondingly smaller signal to noise ra- 
tio). In Fig. 7 we plot the Fourier amplitude of the signal in 
the expected frequency channel (assuming the frequency 
power law with p = 1) as a function of time together with 
errorbars (calculated from 10 neighbouring amplitudes). 
Single data points are clearly very noisy, therefore, not 
much can be said about the time dependence - in par- 
ticular the glitch of July 1996 (Jodrell Bank Crab Pulsar 
Monthly Ephemeris ^) did not leave a clear mark. The 
best linear fit to our data points gives an amplitude of 
Am(t) = (0.0045 ± 0.0006) - (0.00029 ± 0.00026)i, where 
t is expressed in years since 1996.0. 

Does the pulsar free-precess? We presented all the ev- 
idence that we have and are inclined to believe that it 
does point to an affirmative answer. Open questions re- 
main such as: 1) does the amplitude of the free precession 
really change with time? and 2) on what time scale does 
the pulsar relax its internal stress? Observation of pulsar 
free precession can teach us much about pulsar physics. As 
shown by Cadez, Galicic and Calvani 119971 the free pre- 
cession model can quite sensitively (if something is known 
about the state of internal stress) distinguish between dif- 
ferent neutron matter equations of state. Observing the 
change in free precession amplitude and frequency can al- 
low one to learn about external torques on the pulsar (jets 
may produce them) and also about exchange of angular 
momentum between the crust and the neutron superfluid. 
Many questions regarding the 60 g ee period remain open 
and, given the relevance of those questions to our un- 



6 Note that our argument differs from those of D'Alessandro 
and McCulloch 119971 . who rely on Shaham's 119771 super- 
fluid vortex theory in estimating that the angular momentum 



interplay between superfluid interior and the rigid crust is the 
most important mechanism determining the free precession fre- 
quency, and that of Melatos|2D0D] who argues that luj v = euj rot 
"where e is the non-hydrostatic ellipticity" . 
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derstanding of pulsar physics and their interaction with 
neighbouring plasma, it would be important if other re- 
search groups could extend our analysis and observations. 
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